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ABSTRACT 


The DUNE2D boundary layer model for small-scale morphology 
(Andersen, 1999) is compared with bedform morphology measured on the inner shelf in 
11m water depth during the SHOWEX experiment at Duck, N.C. The model consists of 
inter-linked modules for flow, sediment transport and morphology. The flow module is 
based on solving the Reynold’s averaged Navier-Stokes equations in the vertical plane 
with k-omega turbulence closure. The model has been extended to accept a general (but 
periodic) bottom boundary to be able to compare with field data. Boundary layer 
velocity profiles were measured using a Bistatic Coherent Doppler Velocity profiler 
(BCDV). A two-axis scanning sonar altimeter measured small-scale morphology over a 
1 by 1.5 m area with 4 cm horizontal and 0.25 cm vertical resolutions. Bottom maps of 
small-scale morphology were obtained continuously every 20 minutes. A relatively 
simple data sequence was selected for model comparison, during which time the wave 
forcing evolved abruptly from H;j.=0.3m to Hsie=3.0m (bed velocity < 1 m/s), and the 
bed evolved from no motion (relic) to actively migrating vortex ripples. SHOWEX 
bedform changes under low wave plus collinear current conditions resulted in minor 
changes of the yore ripple fields. Bedform migration rates of the model were similar to 
the field migration rates. Like the field data, the modeled data under strong forcing 
removed smaller scale vortex ripples and redistributed the sediment into a larger scale 
ripple with a large portion of sediments in suspension above the bed. Limitations of the 
model owing to the 2-D assumption, periodic boundaries and monochromatic wave 


forcing constraints are discussed. 
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l. INTRODUCTION 


This is a study of ripple migration, evolution, and erosion on the inner continental 
shelf over time. Despite several decades of analytical and empirical research, 
understanding the evolution of bedform morphology remains a fairly young field. In this 
thesis, the morphology problem is considered for non-breaking waves in the shallow-to- 
intermediate water depth range of the inner shelf. Current research emphasis is placed on 
numerical simulation of this problem. This chapter provides motivation for study, 
research objectives, and background on the bottom boundary layer with emphasis on the 
interactions of hydrodynamic forcing on sediment transport and bedforms. 


A. MOTIVATION FOR STUDY 


Naval operations have rapidly transitioned from “blue”, open ocean into coastal, 
shallow and intermediate waters. The operational oceanographer is required to implement 
ocean models of the littoral region to produce accurate forecasts of ocean processes. 


Figure (1.1) lists some Navy operations, as well as dynamic processes that are of concern. 





Figure 1.1. ioral Naval Interests. errr eaiee ti nps. navy. mil/~stanton/miso/) 
] 








Over the past decade, naval research has been highly focused on adapting to this 
new theater of operation. Field studies, laboratory experiments, and numerical modeling 
research of these highly spatial and temporal varying processes are motivation to improve 
physical understanding of the waves, currents, sediment transport, and bottom 
morphology in the littoral region. Advances are slow and expensive. Hopefully, the 
knowledge gained through the combined observations and simulation of the nearshore 


environmental processes will provide accurate oceanographic forecasts for required 


operational applications. 


Common to all dynamic processes in the littoral region are the hydrodynamic 
forces produced through wave, current, and tidal energy. Motivation for this thesis is 
based on understanding wave and current interactions with the bottom boundary layer 
that force sediment transport in the form of bedload and suspended load, and result in 


bedform evolution and migration. 


Wave-formed ripples are the dominant bedforms in coastal regions (Traykovski, 
et al, 1999). Over time, mobile beds can evolve in such a way that the overall bed 
roughness is changed. Furthermore, migration of sediment affects form drag on the 
existing hydrodynamic processes, net sediment transport, and the bottom boundary layer 
turbulence generation. These changes can have significant affect on wave energy 
dissipation (Ardhuin, 2000). 

B. OBJECTIVES 
The objective is to simulate the migration and evolution of small-scale bedforms 


measured in the field using DUNE2D, a numerical model of the bottom boundary layer 


for small-scale morphology. The model is initiated using bedform and velocity 
2 














measurements acquired during the SHOWEX field experiment conducted from 
September to December 1999 at Duck, North Carolina. Small-scale bottom maps and 
high-resolution velocity measurements of a specific storm event are used to qualitatively 
verify DUNE2D output. Efforts are focused on simulating actual field conditions from 
an initial bed configuration. The long-term objective is to provide model validation to 
enhance numerical prediction of nearshore processes, particularly to predict ripple 
formation and migration. 
C. BACKGROUND 

1; Flow Properties 

Accurately determining the flow field in the bottom boundary layer 1s essential to 
predicting sediment transport and morphology evolution. Flow in the ocean is by nature 


turbulent. Sediment transport calculations require estimates of the bed Stress, Tp, 


T, = Pe (1.1) 


where the shear stress is proportional to the frictional velocity, u»’. Flow field dynamics 
are complicated by irregular waves. In the inner shelf regions, flow near the bed is 
influenced by wave groupiness. Quantitative representation is made difficult through 
multiple wave packets interacting with each other superimposed with mean currents of 
varying strength and direction. Deigaard et al (1999) found net sediment transport under 
regular waves and wave groups with bound long waves to be dependent on sediment size, 


as well as wave characteristics. 





2. Sediment Transport through the Shields Parameter 
Sediment transport can be described in terms of the Shields parameter. The 
Shields parameter is the ratio of the destabilizing forces of the bottom stress and the 


stabilizing force of gravity acting on the sediment, 


_ t, 
d= anne “Dd (1.2) 


where g is the gravitational acceleration, d is the sediment diameter, and s is the ratio of 
sediment density to fluid density. Sediments have a critical shear stress that 1s dependent 
upon their grain size, density, porosity, and shape. Once the shear stress exceeds this 
critical threshold, the viscous sub-layer breaks down, and the sand grains are set in 
motion. As the strength of the flow increases, the grains are lifted off the bed. As long as 
the settling velocity of the sediment is small compared with the forces producing lift, the 
grains will remain in suspension above the bed. The grains redistribute back to the bed 
when the fluid forces acting on the grains become less than the gravitational forces. 
Deposition is difficult to predict, as the flow field has eddies of varying strength and size. 
For high Shields numbers (9>1) under non-breaking waves, Diegard et al (1999) found 
that sheet flow occurs, where a large portion of the bedform sediment is set into motion, 
such that the remaining bed is planar in shape. For forcing with low Shields numbers, 
ripple-dominated bedforms occur. 

3. Bedforms 

Ripples are classified by their geometry and commonly divided into two separate 


spatial scales. Dunes, antidunes, megaripples and relic ripples are considered large-scale 


bedforms. Small-scale ripples include vortex, cross-ripples and rolling grain ripples. 














Small-scale ripples are categorized by a length of less than 0.6 meters and a height 
roughly ten percent of the length, or below 0.06 meters (Engelund and Fredsoe, 1982). 

Bagnold (1946) described two types of wave-induced ripples: rolling grain and 
vortex ripples. For increasing Shields numbers, as the critical value is just exceeded for 
flow over a flat, sandy seabed, grains oscillate back and forth in a rolling motion. This 
condition is difficult to find in nature (Andersen, 1999). As the flow intensity increases 
and the turbulent eddies induce more sediment transport, the oscillatory wave action 
organizes the moving grain into the classic vortex ripple triangular shapes. Over time, 
the sediment is moved from the trough of the ripple, and deposited through bedload and 
suspended load towards the mpple peak. The steepness of the mpple increases. 
| Eventually, the flow vortex, formed by the fluid moving over the ripple peak, separates 
from the bed and induces further sediment transport (Neilsen, 1981). Many studies have 
been conducted to find equilibrium or quasi-equilibrium forms of the bed due to wave 
and current forcing with the objective of predicting ripple wavelength, height, and 
steepness as a function of mobility number or Shields parameter. 

4, Numerical Solutions 

Several numerical approaches have been used to describe the flow in the bottom 
boundary layer over an arbitrary bedform. Permier’s (1996) discrete vortex (DV) method 
does a good job describing the qualitative features of the flow. Turbulence models 
provide better quantitative description, which is required for morphological calculations 
(Andersen, 1999). The «x-¢ turbulence closure technique was used by Sato (1986), 


Tsujimoto ef al. (1991), and Perner (1996), while Wilcox (1988, 1993b) introduced the 








K-@ turbulence closure. The x-@ model handles the strong pressure gradients near the 
separation point of the turbulent eddies better than the k-e model. 

The present DUNE2D model study uses the x-a closure. The model is presented 
in the next chapter. Good results have been obtained by Andersen (1999) using 
DUNE2D for monochromatic wave and wave plus current forcing compared with 


laboratory data (Fredsoe ef al., 1999). 


II. DUNE2D MODEL FORMULATION 


DUNEZD is a bottom boundary layer model for sediment transport over small- 
scale ripple morphology (Andersen, 1999). The model resolves two-dimensional flow 
over a one-dimensional bedform. DUNE2D is a modular Fortran code, which evolved 
over many years of research on flow, sediment transport, and most recently, morphology 
at the Technical University of Denmark. DUNEZ2D is composed of three modules: flow, 


transport, and morphology. A flow diagram of the model is represented in Figure (2.1). 





SEX 










4s 
% 
yyy 

SILL SIE MMA nO 

SEE Yi fp LD 
Wj yyy Uy 
LM MY LEE: 
OL YLLLLE LEE 
iE pants mort 


i ithe ef hs 
Ze 
“a 


ARS Noyce. PS See 
AEM AES FAURE RARE ERR ee eee IRS STN TS 
SSRIS SEE ARRAN EARS ASR RR SRR ARCHER CRT REC ECRNSE 
‘! Cae $5 % " fete! aa 
SERINE: SAS REE 







a 
4, 





SSRER SR Sr BEERS 






Sah ES eS I 
SHORES 


SS KP hh 
SSSSEERESA ERS 


5 
3 





Figure 2.1. DUNE2D Model Schematic Flow Diagram. 


The program has several modes of simulation. For each mode, several physical 
assumptions and numerical modeling techniques exist. The model can be run for flow 
and sediment transport simulations, or for bottom boundary layer flow alone. For 
purposes of this study, systematic solutions for all three principle modules are considered 
with emphasis on the evolution of bedforms. Therefore, the model description in this 


chapter will be presented from a morphology simulation perspective. 
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Dynamic forcing is through periodic, sinusoidal wave input. If desired, a mean 
along-axis current can be included. At each grid point of the model, the horizontal, u, 
and vertical, w, velocity flow components are evaluated. Flow solutions are passed to the 
transport module where sediment characteristics, as well as the flow solutions, influence 
sediment transport calculations. Strong turbulent flow is capable of sediment transport, 
both along the bed and in suspension above the bed. Transport solutions are passed to the 


morphology module where bedform changes are integrated over time. 


DUNEZ2D models oscillatory flow with an optional mean collinear current over a 
rippled bedform, as illustrated by Figure (2.2). The twelve parameters listed in Figure 
(2.2) are required to describe the model domain. Appendix (A) provides a further 
explanation of the input parameters. Hydrodynamic forcing, geophysical properties of 
the sediment, numerical schemes, grid definitions, boundary conditions, transport 
methods, and morphology setup are some of the categories that are used to describe 
model input. Appendix (B) provides examples of the ASCII input files and the 


procedural method used to implement DUNE2D code. 


v: kinematic viscosity 
p: fluid density 7 
p,: sediment density -————> 
A: ripple length k 
h: ripple height U, N 
d.): sediment grain diameter i 
D: depth of flow D 
U,: amplitude of oscillatory 
wave velocity 
T: period of oscillatory wave h Ps 
U,: mean current y 
g: gravitational acceleration 
k,: Nikuradse bed roughness C50 


rs 


Figure 2.2. Flow over a nipple bed. 
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A. FLOW MODULE 

The flow module models fully turbulent flows, with Reynolds numbers ranging 
between 10° to 10°, in the bottom boundary layer. The model initializes input variables 
and sets necessary default values, and then solves the flow field given specific dynamic 
forcing over an arbitrary one-dimensional bedform. 

I. Conservation of Momentum and Mass 

The two-dimensional flow equations along the horizontal and vertical axis assume 


an incompressible, Newtonian fluid. The conservation of momentum is, 


” | Ou, 
OB i oye 3 Hs +) i,j=12 @) 





at’ ax, ~—p Ox, ax, |2[( ax, ax, 


where p, uj, and v represent pressure, velocity, density and kinematic viscosit 
> y, 


respectively. The continuity equation is described by, 


—=0 i=1,2 (2.2) 
The Reynolds-averaged Navier-Stokes (RNS) momentum equations are obtained 
by applying Reynolds decomposition and averaging to Equations (2.1) and (2.2) with the 


pressure and flow quantities divided into a mean and fluctuating components (e.g. 
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Ot ON p Ox, Ox,|2\ Ox, Oo, ons 
OU; _ 0 (2.4) 
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where the mean component resolves the waves and currents, denoted by a capital letter, 
the turbulence fluctuations are denoted by primes, and the over bar indicates time 
averaging. A slow time variation of the mean velocity is allowed. 

The last term in Equation (2.3) is the Reynolds stress tensor. The RNS equations 


are closed using turbulent viscosity closure, where the stress tensor is modeled as 


| 2 au, } 
al ae, * Ge 
sade ah Oe ls Shs a ae (2.5) 


ex, oo 


J 


where k= = +yw? ) 


vr is the turbulent viscosity and 6; is the Kronecker delta. Equation (2.5) describes a 





dampening and diffusion process. 


a. Turbulent Kinetic Energy Equations using x-@ Model 


A two-equation K-@ model is used to determine eddy viscosity as a function of 
time and space (Anderson, 1999). The x-@ model (Wilcox (1988; 1993b)) was selected 
over the more commonly used K-E model to accommodate boundary conditions on a 
rough bed in unsteady flow with strong pressure gradients along the bottom boundary. 
The turbulent viscosity term in Equation (2.5) has a direct relationship to « (kinetic 


energy) and @ (vorticity) through a model constant Y. 


np 2.6) 
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The turbulent kinetic energy equations are used to solve for « and @. (Wilcox, 1993b) 





OK OK O O ——<dO ‘ 
U ,-—— = —| (v+ uu.— -BK 
a ax, = Gx, ( ie | Oe. pee (2.7) 
Ow Cw O w — ou, 
| eee (v+ov,) = eal — UU, B : 

Or Ox, ox, ; ox, (2.8) 
where closurecoefficients: y° =1: = Ee B = 2 

ae ED 100’ 

Opts coh tk 
ae Or= >? oO = 7 


These equations are similar to the x-e model, where ¢ is the local dissipation rate 
of turbulent kinetic energy. The term @ does not have a direct physical meaning, but can 


be described as the vorticity of turbulent eddies (Saffmen, 1970), and is related to ¢ by, 


é= pKa (2.9) 


3. Roughness of the Wall 
The bottom boundary is scaled using the concept of a rough wall through the 
Nikuradse roughness, ky/30. In DUNE2D, kn is set to two and a half times the median 


orain diameter. (Fredsoe and Diegaard, 1992). 


ky, = 2.5d,, (2.10) 


The law of the wall is applied to the bottom boundary condition associated with 
the sandy bed. A modification is made to a smooth wall to fit the empirically rough wall 


solution. (Andersen, 1999). The bottom boundary condition is a no-slip condition for the 


1] 





flow, which implies a viscous sub-layer and x=0 at the bed, while @ follows (Wilcox, 





1993b): 
2 
Vv 
where U, = friction velcity = jaa 
fe (2.12) 


— [(50/Kg)?, kez <25 
- | ae kez >25 
Andersen states that the methodology for attaining a rough wall solution through 
parameterization constant, Sr, is not physically correct. Wall functions are normally used 
with turbulent kinetic energy models that don’t involve viscous sub-layers. The x-o 
model has a viscous sub-layer due the no-slip condition imposed at the bottom boundary. 
Andersen also points out that model results using the k-o model are sensitive to the free- 
stream boundary condition for @. To overcome this problem, a symmetric boundary is 
applied to the upper boundary. (Andersen, 1999) 


B. SEDIMENT TRANSPORT 


DUNE2D separates sediment transport into bedload and suspended load 
calculations using techniques described by Engelund and Fredsoe (1976) and Fredsoe and 
Diegaard (1992). Turbulent flow causes destabilizing forces upon the rippled bed. 
Sediment is thrown into motion when strong destabilizing forces, specifically, drag and 


lift, overcome gravitational stabilizing forces. Segregation of bedload and suspended 


load is imposed by an imaginary boundary. Moving sediment in repeated contact with 
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the bed is known as bedioad transport. Suspended load transport is sediment in motion 


above the bedload transport boundary, defined at a height 2d, by a function 


Cy (6, /) z=2d os ra 


where @ is the Shields parameter and y is the slope of the bed and d is the median grain 


diameter. Both transport load calculations are based on the Shields parameter, 0, for a 


flat bed (Equation 1.2), where the maximum shear stress during a wave cycle is 


ml 
t, =— pf, 
an) A (2.14) 
where U, is the wave velocity amplitude and the wave frictional coefficient 1s, 
A —0.25 A 
| ood A > — > 50 * 4245} 
k N N 
where A is the orbital amplitude of the wave motion, 
Ane (2.16) 


20 


The local bed shear must be greater than the critical bed shear for sediment motion to be 
initiated. Bedload is a function of the bed shear, while suspended load is a function of 
bed shear as well as settling velocity. The settling velocity of the sediment is determined 


using an empirical formula for the drag coefficient of sand (Fredsoe and Deigaard, 1992): 


fi (2.17) 
con 
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where Cp is the drag coefficient and Reg is the sediment Reynolds number. 

1. Bedload Calculations 

Bedload calculations are based on the Engelund and Fredsoe (1976) formulation. 
Sediment velocity is found from the balance of the gravity, drag, and friction forces on a 


single grain on a sloping bed. The non-dimensional single grain velocity is given by, 


cy 


9 a =10 (emperical constant), (2.18) 





th wal 1-07 
U, 








where Uy is the grain velocity and 6, is the critical Shields number for a sloping bed, 





sin 7 
6. =0,) cosy + 
: : y ar (2.19) 


D 


in which and o is the sign of the direction of the bed load and the critical Shields 


parameter for a flat bed 1s, 


©  3C,a° 
for which up, is the dynamic friction coefficient. (Andersen, 1999) 
It is assumed only a single layer of grains can be in motion at the bed. From this 
assumption, an equation, based on empirical data, is used to solve the fraction of grains 


that are in motion per unit area, 


7L: 
1 614 (2.20) 














Knowing the number, or probability of grains in motion, bedload transport, or bedload 
flux, the rate of transport of a volume of sediment per time, is determined. The non- 


dimensional bedload flux is (Fredsoe and Diegaard, 1992) 
I 
¢,(0.7)= n—@U, vs — Vga° (2.21) 


After substituting for n, Up , and including Fredsoe’s linear gravity correction to the 


Meyer-Peter flat bed equation for bedload flux, Equation (2.21) simplifies to: 


$,(0,7)= 80-9, - Br)” none 
where § = 0.1. 

23 Suspended Load Calculations 

Suspended sediment is advected by wave motion and mean currents while falling 
towards the bed based on the sediment’s settling velocity, ws. Strong turbulent 
fluctuations with small settling velocity sediment produce sediment in suspension above 
the bedload boundary. Suspended sediment concentration, c, is modeled by, 

c= fed] . ¢, =0.65gm~ (2.23) 


’ 
where 2g, the mean free path of the sediment at the bedload boundary condition, is 


“4 
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Ay 





Ak TE: 3 
4, =——— 9-6, - n=’ : 
— a. ae io , Gay) 


and n, again, represents the fraction of grains in motion per unit area as described by 


Equation (2.20). 
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The total instantaneous suspended concentration is found using Equation (2.23) 
where the time rate change of concentration is calculated using a transport-diffusion 
equation: 


DE any“ (ENe) (2.25) 
Dt Oz 


where the settling velocity is represented by ws using Equation (2.17) and the term €; vr 
is the suspended sediment diffusivity, which is equal to the eddy viscosity defined by 
Equation (2.6). 

The instantaneous suspended sediment flux is obtained through the instantaneous 


velocity distribution and concentration found from Equation (2.25), 


q, = u(x, z,t)e(x, 2,2) (2.26) 


Finally, The total sediment load flux, q:, through a cross section, is the vertical integral of 
the flux profile, 


g.(xth=9,(x%,t)+ [9,(% 2,02 (2.27) 


h(x,t} 


Andersen states that one of the fundamental problems with the calculation of the bedload 
transport is that it is based on mean bed shear stresses. The model does not incorporate 
the turbulent fluctuations, x. Despite this deficiency, dynamic ripple mechanisms are 


preserved. (Andersen, 1999) 
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Cc. MORPHOLOGY MODULE 


The morphology module is based on the continuity of mass for sediment. From 
the total transport Equation (2.27), the change of the bedform over time is found by time 
integrating, 


GH(x,t) __ 1 Oa) (2.28) 


Ot l-p & 
where p is the porosity of the sand. The main interest for the morphology simulations is 
the mean bed height at the end of each wave period. For further discussion of the output 


parameters see Appendix (A). 


D. WAVE PLUS CURRENTS SUPPLEMENT 


The model forcing can be supplemented with the addition of a collinear current on 
the oscillatory flow. The addition of the current is made through a replacement of the 
constant friction factor for wave-only flows to a wave plus current friction factor, fie, 


| U.DY 
ie, at oa (2.29) 
fe : 


where fe is the constant friction factor for the current in a wave plus current flow. Dis 
the depth scaling term and Q is the water discharge per unit width (Fredsoe et al., 1999). 


An equivalent Nikuradse roughness, Arye, is related to fuc by 





Pe 148Dexg =P (2.30) 


V Src 


The advection of sediment by the mean flow is calculated in the bedload and suspended 


load modules. 
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II. FIELD DATA 


A long-term inner shelf study of the bottom boundary layer was part of the 
Shoaling Waves Experiment, SHOWEX, at Duck, North Carolina (Figure 3.1). Multi- 
scale measurements of the hydrodynamic interactions over sandy rippled bedforms were 
made. This chapter describes the data collection procedure, measurement techniques, and 


bedform morphology due to ocean forcing caused by a specific storm event. 
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Figure 3.1. U.S. Army Corps of Engineers Field Research Facility at Duck, North 
Carolina. Site of SHOWEX 99. (from http://www. frf.usace.army.mil/) 





A. MEASUREMENTS 
SHOWEX was an Office of Naval Research funded Directed Research Initiative 


(DRI) to improve wave propagation models in the coastal regions. One of the objectives 
was to measure the dissipation of shoaling waves across the inner shelf with emphasis on 


the affects of the feedback between wave induced ripples and the resulting increased 
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dissipation of the waves at the bed. Measurements of the bottom boundary layer were 
made at 11 meters water depth from a bottom mounted platform with precision 


instruments at site S-3 (Figure 3.2) at Duck, North Carolina, during 28 September to 10 


December 1999. 
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Figure 3.2. SHOWEX 99 measurement platform location along Duck, N.C. (after 
http://cheyenne.rsmas.miami.edu/duck99/ ) The platform is approximately 1.4km off the 


beach. 


Sediment collection of the relatively homogeneous fine sand was made in the 





observation area four times during the three-month experiment. The median sediment 


diameter was 0.1 mm based on a sieving procedure. 


A schematic representation of the instrument platform is shown in Figure (3.3), 
illustrating the sensors used to measure velocities within the bottom boundary layer as 


well as the bedform and suspended sediment concentrations. A 1.3 MHz pulsed, bistatic 
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Figure 3.3. Bottom boundary layer measurement platform. 


coherent Doppler velocity (BCDV) profiler (Stanton 1996, 2000) measured vertical 
profiles of three component velocities, u, v, and w, and sediment concentration in 2 cm 
diameter volumes of water every 0.6 cm over the bottom 60 cm, starting at 30 cm from 
the sensor face. The sampling frequency was 18Hz, resulting in high resolution velocity 
profiles. The BCDV’s orientation is illustrated in Figure (3.4). Tilt sensors in the BCDV 
allowed precise orientation of the vertical coordinate system. The raw BCDV vector 
velocity profile data were transformed into an earth-referenced coordinate system during 


post-processing. 
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F:gure 3.4. Platform measurement orientation is 22° to the left of Beach North. 
Piatform located in 11m depth. 


A two-axis, scanning, high-resolution acoustic altimeter was collocated with the 
BCDV profiler on the instrument platform. The scanning altimeter mapped a 4m by 2m 
bed area every 15 to 30 minutes. Resolution of the bed was 4cm in the horizontal and 
0.25 cm in the vertical over an inner 1 m? area. Horizontal resolution decreases away 
from the measurement center due to step sensor angles and distances of travel. Post- 
processing of the raw range data included removal of false targets, transformation of 
range and tilt angles into Cartesian coordinates, and objective analysis techniques to 
smoothly map the scanned data points across the measurement area. Final bedform 
observation boundary limits, shown in Figure (3.5), were reduced to 1.5m by 1m box 


defining an area of maximum resolution. 
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Figure 3.5. | A two-dimensional bedform mapped using X-Y acoustic altimeter during 


SHOWEX ‘99. Black line is the BCDV measurement center. Right foreground is a hole 
generated by a live shellfish embedded in the ocean floor. 


B. VELOCITY FORCING AND BEDFORM MORPHOLOGY DURING 11-13 


NOVEMBER 1999 


The evolution and migration of bedform are examined for a storm event during 
year days 315.5 to 317 (November 11-13). The storm produced strong local winds that 
caused the significant wave height to increase from 0.3m around year day 315.5 to a peak 
of nearly 3m at year day 316.4, as measured by the FRF 8m bottom pressure gauge array 
(Figure 3.6). During the storm, wave forcing increased the orbital displacement just 
above the bottom boundary layer from 0.1m to 1.1m. As the storm tracked out of the 
region, local wind waves and swell dissipated. 
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Figure 3.6. | FRF Duck, N. c measured weather and oceanographic data from archived 
data collection for November 1999. Box indicates the storm event of 10-13 November. 
(after http://www.frf.usace.army.mil) 


The uw and v velocity components measured by the BCDV were arithmetically 
averaged over the vertical between 10cm to 30cm above the bed to form 278 second time 
series at 0.1 year day intervals (2.4 hours). The velocities represent the near-bed 
velocities just above the turbulent boundary ave: sometimes referred to as ua and Vo. 
From the wu and v time series, the wave oscillatory velocity amplitude, Uo, peak wave 


period and mean currents were determined (Figure 3.7). 


24 








—— U. (Wave Velocity Amplitude) | gomnnnnngna 
se= U (Along-shore Current) 2 , 
{| =" V_(Cross-shore Current) 
9.67) wom 7 (Peak Wave Period) 


S 
= 


o 
N 
Peak Wave Period (s) 


Hydrodynamic Velocity (m/s) 


Rl 
rd 


Ot met my, 
of 


: tia . zs" 

0 ‘to nga 
# 

é 

F oe” : o” 
- ee 

ae gn ag 

ag ai a 
7 ig i el ee 





oy ieee | | 0 
315.4 315.6 315:8 316 316.2 316.4 316.6 
Year day 
Figure 3.7. Meancurrents and oscillatory wave velocity amplitude for year days 315.5 
to 316.6 illustrating storm event spin-up of the ocean dynamics. 


The wave velocity amplitude is calculated.as the equivalent sinusoidal velocity 


amplitude that gives the same wave variance as measured, 


U,=2Jo2+0 (3.1) 


where 6,” and o,” are the along-shore and cross-shore velocity component variances. 
The peak wave period of the initially mild waves, composed of local and distant swell 
was not well defined. By year day 315.9 the storm peak wave period became well 
defined and increased over the following 24 hours. The majority of the mean current 


energy was in the along-shore direction as the magnitude increased from 15 cm/s to a 
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scanned acoustic altimeter data are shown in Figure (3.8). Initially at year day 315.5, the 
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Evolution of bedform along centerline, cross-shore transect during storm 
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IV. DUNE2D MORPHOLOGY OUTPUT ANALYSIS 


Sensi tests of DUNE2D morphology simulations were first performed to 
determine model output tendencies given specific monochromatic wave forcing and wave 
plus current forcing. The model morphology evolution and migration is then compared 
with field data observed during a storm event of SHOWEX ‘99. In both studies, a cross- 
shore transect bedform was chosen from the SHOWEX ’99 altimeter data and 
transformed into a periodic boundary bedform represented in Figure (4.1). All 


simulations use fine sediment grain diameter equal to 0.1mm, as measured at SHOWEX. 





0 02 o4 06 £08 1 12 14. 16 £18 
x (m) 


Figure 4.1. | Transformed cross-shore transect of bedform from SHOWEX °99 data 
(year day 315.9) used to initialize DUNE2D simulations. The beach is to the right with 
the initial oscillatory wave forcing going towards the right (See Appendix (C) for 
explanation of bedform transformation). 


A. MODEL SENSITIVITY TESTS 


Sensitivity studies of DUNE2D morphology simulations are conducted on the 
initial bedform (Figure 4.1) using varying initial wave forcing and wave plus current 
collinear forcing over a short simulation time duration of 20 wave periods to obtain a 


qualitative feeling for the model performance. 


pe | 





1, Wave-only Case Study 


Four wave-only forcing tests were made with different monochromatic wave 
input values as summarized by Table (4.1). Tests 1 and 2 represent the lower wave 


Table 4.1. Four test cases in wave-only sensitivity study of morphology simulations. 


Oscillatory — Uo 0.35 0.35 0.70 0.70 


Orbital Displacement A, (m 


Mobility Number, wv 


76 76 303 303 
Field Measured Bedform 315.9 315.9 315.9 
(Year day) 


energy simulations, with a mobility number, w, equal to 76, indicative of the pre-storm 
condition. Tests 3 and 4 are the higher wave energy simulations occurring during the 


storm event, with a y=303. 


Bedform changes over twenty wave periods for varying monochromatic 
oscillatory wave-only forcing is shown in Figures (4.2a-d). In each figure, the output 
bedforms for the first and twentieth period are displayed. Similar bedform changes 


resulted for test cases 1 and 2, while test cases 3 and 4 are similar. Vortex ripples are 





created over top the larger scale bedform in test cases 1 and 2. The amplitude of the 
vortex ripples is greater for the longer period waves of case 2. The wavelengths are the 
same for both cases with the average distance between the new ripple peaks is equal to 


10cm. 
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Test #3. 
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Figure 4.2. | DUNE2D bedform morphology for field bedform year day 315.9 with 
wave-only monochromatic forcing after 20 wave periods. (a)Test#1: Up=0.35m/s, 
T=4.79s, w=76; (b) Test#2: Uo=0.35 m/s, T=10s, y=76; (c) Test#3: Up=0.70m/s, 
T=4.79s, y=303; (d) Test#4: Up=0.70m/s, T=10s, yw=303. 
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For the high mobility number cases 3 and 4, the vortex ripples are planed off, 
leaving a smoother large-scale bedform. The model tends to quickly change the bed after 


a short time equal to 5-10 wave periods. 


In summary, low energy wave forcing used in tests 1 and 2 results in the creation 
of vortex ripples on top of the larger scale bedform. The high energy wave forcing used 
in tests 3 and 4 planes off small-scale ripples on top of the larger scale ripple bedform. 


Zz Wave plus Collinear Current Case Study 


The effects of adding a collinear current to the same low amplitude, short period 
wave forcing and initial field measured bedform (SHOWEX ’99 year day 315.9) are 
examined. Four separate collinear currents (Table 4.2) are considered where the ratio of 
current to wave velocity amplitude is U/U, = 0.5 (Test 5), U/U, = 1.0 (Test 6), U/U. = 
2.0 (Test 7), and a reverse direction current, U,/U, = -2.0 (Test 8). Again the model was 


Table 4.2. Model hydrodynamic input for wave plus current sensitivity study of 
morphology simulations. 
n 


Py lw Tt | lt 


Field Measure Bedform 315.9 315.9 315.9 315.9 
Year day 


run for twenty wave periods for each test case with the results shown in Figure (4.3a-d). 












The results of test case 1 for the wave only forcing are included in the figures for 


comparison. 
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Test #5 





Figure 4.3. | DUNE2D bedform morphology given an initial bedform (field year day 
315.9) with monochromatic wave (Up=0.35m/s, T=4.8s, y=76 ) plus current forcing after 
20 wave periods. (a)Test#5:U/U.=0.5; (b) Test#6 U/U.=1.0; (c) Test#7: U/U.=2.0; (d) 
Test#8: U/U,=-2.0. 
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For U./U, = 0.5 (Figure 4.3a), the current has little effect on the ripple bed. The 
20" period bedform for test #5 is very similar to the bedform morphology of the wave- 
only test case, which supports Andersen’s (1999) model results for wave plus low current 
magnitudes (U./U, < 0.5) having minor effect on bed morphology. For U./U. = 1.0 
(Figure 4.3b), there is a noticeable migration of the bed in the direction of the current. 
For U,/U, = 2.0 (Figure 4.3c), more sediment is removed from the bed and placed into 
suspension and greater migration of the bed. Finally, reversing the current flow direction 
180 degrees, shifts the bed migration to the opposite direction. In summary, a collinear 
current's effect on the bed is significant at values where the current magnitude is equal to 
or greater than the wave velocity amplitude, U./U, >1.0. 


B. MODEL COMPARISON WITH SHOWEX ’99 OBSERVATIONS 


A DUNE2D morphology simulation of the SHOWEX 99 storm event described 
in Chapter III was conducted to assess how well the model compares with real-world 
bedform evolution and migration. In a course of the day, field-measured bedforms 
experienced low wave energy initially, followed by large wave energy generated by a 
coastal storm system that resulted in a modified bottom bedform. The bedform evolved 
from a mega-ripple with superimposed small-scale roughness created by vortex ripples to 
a flat, smooth bed twenty-four hours later. Refer to Figure (3.7) for the field bedform 
morphology time series produced by this storm-forcing event. . 


1. Simulation of Storm Event 


The storm event is simulated by two separate model runs representing the initial 


pre-storm, low wave and current energy, followed by the conditions at the peak wave 


energy plus current during the storm. The onshore, cross-shore current was low for both 
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simulations. Most of the current was oriented in the along-shore direction. The 
orthogonal currents can not be directly included into the simulation since DUNE2D is 


limited to collinear currents. 


An approach to incorporate the effects of the along-shore currents on the bedform 
is to use an effective wave velocity amplitude that is equivalent to the combined wave 


and orthogonal current bed shear stress, 


(4.1) 


1 effective = ty = ie 


where the wave bed shear stress is calculated using Equation (2.14) and the orthogonal 
current stress is calculated using Soulsby’s (1997) equation with a current only drag 


coefficient. 


(4.2) 


From the field measurements, the calculated along-shore, current-only shear stress is 1/3 
the wave-only bed stress. Therefore, an effective stress of 133% of the wave-only stress 
is used to include the effects of the orthogonal current. An effective wave velocity 
amplitude, U., is obtained from this effective bed stress using Equation (2.14). The 
results of the modeled, morphology bedform evolution using this effective wave bed 
shear stress failed because the model generated improper equilibrium ripple lengths and 


profiles. Therefore, the along-shore current effects on the bedform are assumed not to be 
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important, which is supported by laboratory results that indicate orthogonal currents do 
not significantly affect maximum shear stresses under oscillatory flow (Soulsby ef al, 
1993). 

Long-time simulations were made to ensure the bed would reach equilibrium, 
such that accurate ripple geometry measurements could be compared to the field 
bedforms. The lower energy simulation (Table 4.3., Test #9, y=16) represents the wave 
conditions at the beginning of the storm, near year day 315.9. | The high-energy waves 
conditions on year day 316.4 (Table 4.3., Test#10, w=341) are used to force the initial 
bedform of year day 315.9 to examine the morphology evolution and migration during 


the peak of the storm. 






Table 4.3. Two test cases in storm event case study of morphology simulations. 


Hydrodynamic Input Test #9 (Low Enereg Test#10 (High Energy) 


Uo (m/s 






P 4} 
=} Fa) 






Uc (m/s 


Bedform (Year day) 


Total Simulation Time (s) 





The simulated bedform morphology time series are shown in Figures (4.4) and 





(4.5a) for the low hydrodynamic energy with the field bedform morphology evolution 
following the same forcing given in Figure (4.5b) for comparison. The modeled 
morphology for the high energy simulation (test #10) is illustrated in Figures (4.6) and 


(4.7a), with the field bedform results in Figure (4.7b) for comparison. 
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Figure 4.4. Test#9: DUNE2D morphology time series given an initial bedform (year 
day 315.9) with low wave energy (Uo=0.16, T=4.8, Uc=-.11, y=16, A=0.12m). 


Period # 1 (+0.00 hrs) 
«+ Period # 200 (+0.26 hrs). 

Period # 1200 (+1:60 hrs) | / 

Period # 2000.(+2.66 hrs) #7 /' 





Figure 4.5a. Test #9: DUNE2D morphology bedform snapshots over simulation time 
resulting in a vortex nipple field at equilibrium state (ripple geometry: Aave=0.11m, 1 ave 
(height) = 0.02m, 1 avg /A ave (Slope) = 0.18, Aave/A= 0.92). Migration: 4.2cm/hr onshore. 












—— Year Day:315:8979 (4+0.00:hrs) 
~~ Year Day:315.9344 (+0:88 hrs) 
~~ Year Day:315.9553 (+1.38 hrs) 
ey | --- Year Day:315.9991 (42.43 hrs) 









~.. ow gh A a en 
Meet “ o 


~~ 


Bed Migration 


0.5 


x (m) 
Figure 4.5b. SHOWEX ’99 Bedform morphology time series snap shots taken from 
year day 315.9 to 2.4 hours later. (superimposed vortex ripple geometry: Aave=0.2m, 
Nave = 0.008m, 1 ave /A ave = 0.04) Migration: 4.5cm/hr onshore 


35 





SRO 


Os 


‘ 
ter 
oe *8ay 
tony 
ay 


eee 6 be bee ag weet 
oreenee 
peat enerr® 


; : “S 
“Oty Eee S ae aie wee ate hy 
: ely sees eee ee ie cern cemnren nes’ THe, 
teea “8S ene arty ala a a hs, 

: oe ; 

! sean, 4 , : 
oF es . Pe... . 
7 es dee, ‘ i 

cee tO : : ‘ 

ree, vey ‘ 
: wets 


ee, 
; 
teen, 
e ‘ + 
rr thea ego ances 









toot sae 4 
ate FOSS SSANCEC ence 


Period # 


x.(m) 
Figure 4.6. | Test#10: DUNE2D morphology time series given an initial bedform (year 
day 315.9) with high wave energy (Uo=0.74, T=9.6, Uc=-0.1, y=341, A=1.1m). 


0.06 
— Period # 1 (+0.00 hrs) 
~<«- Period # 200 .(+0.53 hrs) 

Eo. 04] —<— Period # 600 (41.59 hrs) 
--= Period # 1000 (+2.66 hrs) 


Figure 4.7a. Test #10: DUNE2D Model morphology bedform snapshots over 
simulation time resulting in a flat mega-ripple equilibrium state (ripple geometry: 
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Figure 4.7b. SHOWEX ’99 Bedform morphology time series snap shots taken from 
year day 316.4(bedform during peak wave forcing) to 3.6 hours later. (Ripple geometry: 
Aaverage=1.5m, 1 = 0.02m, n/A = 0.02). Migration: 2.2 cm/hr on-shore. 

36 





2: Discussion of Storm Event Simulation 


The low energy simulation (test #9) reached equilibrium after about 2000 wave 
periods (+2.6hrs). The bedform is characterized by small-scaled vortex ripples over top a 
larger scale ripple, much like the field bedform (Figure 4.5b). The average wavelength of 
the vortex ripples at equilibrium state is 0.11 m, with an average height of 0.02m. A ratio 
of average ripple wavelength to orbital displacement is approximately equal to one. This 
compares qualitatively well with the observed bedform during the low wave energy state 
and with Neilsen’s (1981) equation (3.4.3) for V2(A/ A), given w<200. The slope of the 
measured vortex ripples superimposed over the larger scale ripple (Figure 4.5b) is 0.04 
compared with the modeled ripples slope of 0.18. Some of the difference is due to the 
field instrument’s 4cm spatial roll-off in the horizontal that flattens off the top of the 
small-scale ripples of the bed. However, since the observed forcing had a finite spectral 
width and directional spread, the SHOWEX forcing does not closely match the single 
frequency unidirectional model forcing. The migration of the modeled vortex nipple is 
4.2 cm/hr onshore, while the field bedform is similar at 4.5 cm/hr onshore. A summary 
of the results is given in Table (4.4). | 


Table 4.4. Storm event ripple geometry comparison (DUNE2D vs. SHOWEX 
Ripple Geometry 
Output DUNE2D SHOWEX DUNE2D SHOWEX 
Test #9 (YR 315.9) Test #10 (YR 316.57) 
= 13 


0.008 | _0.008 
Slope (naval Rave) 0.004 


hg /A(t ee 14 


Migration (cm/hy) ee ee 












The high-energy simulation (test #10) reaches equilibrium in about 1000 wave 
periods (+2.6 hours), faster than the field bedform (+3.6 hours). The average wavelength 
of the flat, large-scale ripple is 1.8 m, with a height of 0.008 m. The ratio of ripple 
wavelength to orbital displacement is around 1.5. The simulated large-scale, flat ripple 
has a similar dip in the middle of the bed with higher peaks towards the ends of the 
boundaries as the field bedform (Figure 4.7b). The modeled migration rate is unresolved, 


while the measured bed, under the high energy forcing, migrates at 2.2 cm/hr onshore. 


A question arises with the ripple geometry for the high energy simulation (Test 
#10) as to whether its equilibrium bedform is influenced by the domain size of the grid, 
as the computed equilibrium ripple wavelength, Aave, appears to be equal to the 
horizontal domain of the model run. (Figure 4.7a) To examine the question, four separate 
model runs were performed where the length of the domain of the initial bedform was + 
varied, as summarized in Table (4.5) and Figure (4.8a-d). 


Table 4.5. Grid test cases to verify boundary effects on morphology simulations. 


Hydrodynamic Input Test #11 Test#12 Test#13 Test #14 
Uo (m/s), T(s), w 0.69 , 96 , 292 


Initial Bedform (Year day) 315.9 | 
DUNE2D Bedform Horizontal ‘s ‘4 ‘ ; 
Distance Increase from Original Bed L007, 3, SVU 
Total Simulation Time (s) 
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Figure 4.8. | DUNE2D initial bedform for (a)Test #11,(b) Test #12, (c) Test #13 (d) 
Test #14. 
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The high energy model simulations, with similar forcing as test #10, produce 


different equilibrium ripple wavelengths than measured (Figure 4.9a-d). The model 


Yeast €11 


Test #12 


‘Antal, Fe, 





O:CG8: pwr rr itiaad Besct 3 
i awe JoOGoth Poriod Bed i Test#13_. 





- ——7 


. $e . - se a . ; S— ‘ 


Figure 4.9. _ DUNE2D morphology simulations illustrating bedform evolution at 
equilibrium (+2.6hrs) for (a)Test #11(Aave=1.8m), (b) Test #12 (Aave=2.2m), (c) Test #13 
(Aave=1.5m). (d) Test#14 (unresolved Aavzg ) 


output for the 30% and 50% bedform enlargement produce single mode bedforms with 
wavelengths equal to the total horizontal domain, 1.8m and 2.2m, respectively. A 
doubled domain with two bedforms produces two flatten ripples — an average 
wavelength equal to 1.5m, about half of the total horizontal domain. Test #4 has an 
unresolved equilibrium wavelength due to the large perturbation at the domain ends. The 


four tests show the equilibrium wavelength of the large-scale ripple that evolves under 
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high energy forcing (forcing capable of removing vortex ripples) is affected by the 
horizontal domain of the initial bedform. Therefore, the high energy results can be only 
qualitatively compared with the field bedform. The modeled equilibrium bedform show 
similar removal of vortex ripples and a flattened, large-scale ripple for high energy 
forcing as the field measurements. 


c. WAVE SIMULATION OF FIELD STORM EVENT 


DUNEZD is limited to monochromatic wave forcing. However, the field data 
show narrow-banded wave forcing with some groupiness present. To simulate the 
eroupiness of the field measured waves, a sequence of monochromatic wave model runs 
was performed. Each successive model run was represented by the equivalent wave 
velocity amplitude and period of the measured waves. The duration of each run was 
equal to the measured duration of the section of real wave condition, so that a piecewise 
narrow band forcing time series could be constructed. Additionally; each model run was 


initiated by the previous final morphology bedform output and flow solutions. 


However, the method failed when trying to initialize the successive runs with the 
previous bedform output and flow solutions. For the case of a moveable bed as used in 
these simulations, the initial bedform is normalized by a scaling parameter, D, which is 
based on the monochromatic forcing input (U, and T). The model generates a grid from 
this initial, normalized bedform. When the successive model runs are initiated with the 
previous output bedform and flow solutions, the model is unable to resolve the new flow 
solutions due to different scaling values used to generate dissimilar grids between the two 


model runs. Further work is required to incorporate the amplitude changes associated 


with groupiness waves that are representative of the real world. 
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V. CONCLUSION 


The migration and evolution of small-scale bedforms was simulated using 
DUNEZ2D based on field hydrodynamic forcing of the bottom boundary layer. The model 
simulates flow, sediment transport, al morphology in the bottom boundary layer to 
enable bedform comparisons with the field data. The initial cross-shore bedform used in 
the model closely matched the measured SHOWEX acoustic altimeter data. However, it 
is required to manipulate the measured bedform to conform to the cyclic boundary 
requirement. At low energy forcing, the model generated vortex ripples. The 
equilibrium wavelengths of the rippled bedforms were near the orbital diameter for the 
oscillatory wave forcing. SHOWEX bedform changes under low wave plus collinear 
current conditions resulted in minor changes of the vortex ripple fields. Bedform 


migration rates of the model were similar to the field migration rates. 


DUNEZ2D simulations of high energy forcing over a ripple bed resulted in greater 
change in the bed over a shorter duration of time compared with the lower energy forcing 
simulation. Like the field data, the modeled data under strong forcing removed smaller 
scale vortex ripples and redistributed the sediment into a larger scale ripple with a large 
portion of sediments in suspension above the bed. Final equilibrium bedforms consisted 
of larger scale, flattened ripples with wavelengths reaching beyond orbital excursion 
values. However, the horizontal domain of the initial bedform and the number of initial 
low-mode ripples was found to be somewhat influenced by the equilibrium ripple 


wavelength of the large-scale ripple that evolved under high energy forcing. 
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In nature, sediment sizes and distribution varies spatially and temporally. The 
model simulations are limited to a single sediment grain diameter designation, which in 
tum affects the Shield parameter calculations. Ultimately, the Shield’s parameter 


calculations relate to sediment transport and bedform morphology accuracy. 


The oceanic wave forcing is not monochromatic and unidirectional, so detailed 
comparisons with DUNE2D model are not possible. However, the modeled bedform data 
was compared qualitatively with changes in the measured bedforms (Table 4.4). 
Differences can be ascribed to model input limitations, which include simple sinusoidal 
wave forcing, collinear current directionality, periodic bedform boundaries, and single 
sediment grain diameter selection. Efforts focused on generalizing the wave forcing to 
simulate actual field conditions and initial bed configuration. Simulation of the storm 
event through sinusoidal flow over bedforms restricted actual simulation of real waves. 
The higher energy SHOWEX velocity fields were characterized by groupy waves. Wave 
simulation could not duplicate the groupiness observed in the field. The inability to vary 
current direction within the model reduced proper sediment transport. During the storm 
event, strong along-shore current existed. Simulation of collinear currents demonstrated 
substantial transport when currents are greater than the velocity amplitude, U./U, > 1.0. 
Despite the inability to establish arbitrary wave forcing and directional currents, the 
model bedform tendencies based on the simple forcing provided good qualitative 
asreement with field bedforms. Continued work on the capability of the model to 
incorporate real world wave and current forcing is encouraged. The model provides a 


good basis for which to better understand the probable bed migration and evolution over 


time given near monochromatic wave forcing in the nearshore. 
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APPENDIX A. DUNE2D MODEL SOLUTION 


This appendix provides detailed discussion of DUNE2D model solutions. The 
following sections are divided between model input and output descriptions. 


A. MODEL INPUT 


Hydrodynamic forcing, geophysical properties of the sediment, numerical 
schemes, grid definitions, boundary conditions, transport methods, and morphology setup 
are discussed below. 


1. Hydrodynamic and Geophysical Input 


The model is forced by monochromatic sinusoidal wave motion. Wave motion is 
described by the amplitude of the oscillatory wave velocity, Up», and period, T. 
Implementation of a current in the same direction as the wave is an option. The scaling 
depth; D, is applied to most input terms to provide non-dimensionality where 

D = 2m =U,T eel) 


The model is based on turbulent boundary flow. The flow is characterized by the 


Reynolds number computation as described below: 


nee. (A.2) 


e 3 


V 


v = kinematic viscosity 

Sediment description is limited to a homogeneous bottom type, described by the 
medium grain diameter, dso. The grain diameter is incorporated into the model through 
the roughness of the bottom, k, (Equation 2.10). (Fredsoe and Deigaard, 1992) Porosity 
of the sediment is set at 0.4. The ratio of sediment density, ps, to water density, p, 1s 


constant at 2.65, based on in-situ sand samples. 
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Model simulation time is based on total number of wave periods. Time steps, At, 
are a function of the number of iterations for each wave period. Both time functions are 
described by the DUNE2D model input variable WAVE. 


Table A.l. | Hydrodynamic and geophysical input. 


DUNE2D Physical Parameter 
Variable 


ROUG Non-dimensional Nikuradse constant, k,/D 


ANGU Non-dimensional depth, D/a 


WAVE Total number of wave periods & iterations per period 
Reynolds mumber, Re 


2. Numerics and Scheme Setup 












A finite volume numerical method is used to discritize the general flow equations 
in section (2.3). The flow is assumed to be incompressible. The pressure term in the 
Reynolds-averaged Navier-Stokes (RN S) equations is calculated using the PISO 
algorithm (Patankar, 1980). The PISO algorithm is an iterative process used to determine . 
the pressure values at each grid point. The user provides the number of times the PISO 
algorithm is implemented. For flows that are difficult to converge, a higher number of 
iterations are required. Andersen emphasizes the difficulty of solving the pressure term 


within the RNS equations. 


Zijlema’s (1996) ISNAS scheme, a high order spatial discretization method, is 
used to remove numerical diffusion errors. Implicit time discretization methods are 


utilized, except for the advection terms where semi-implicit techniques are implemented. 
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Table A.2. | Numerical Input. 


Physical Parameter 
Variable 


| Discretization using ISNAS scheme _ 
PISO PISO algorithm parameter 
K-@ turbulent closure model 










3. Grid Input 


For morphology simulations, the grid is generated using a transfinite interpolation 
technique. The transfinite interpolation technique produces finer resolution at the bed, 
while vertically stretching the grid points out as the distance from the bed increases. A 
specific constant, GMET, is used to define the degree of stretching. Large GMET values 
correspond to increased stretching. Additionally, the total number of vertical coordinate 


points, M; is a function of the transfinite interpolation method. 


The user can specify bedform through an ASCTI file or have DUNE2D generate a 
bottom using several predefined bedform shapes. For user-specified bedforms, the 
ASCII file contains information about the number of horizontal, N, and vertical, M, grid 
points. Also, the file contains the bedform grid values normalized by the scaling term, D. 


An example of a grid file, test.grd, is given in Appendix (B). 


The resolution of the grid spacing is set by the bed roughness, which is a function 
of the bedform profile, Nikuradse constant, sediment grain diameter and wave forcing 
input. From Fredsoe and Deigaard, the relationship used to determine the resolution 


parameter, RESO, requires the calculation of the maximum frictional velocity, u, 
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k 0.125 
u= 0.1/2U, —2 
a 

(A.3) 


ic 


RESO = integer 227. es 
04 " v 


RESO is used in the grid generating procedure. Large RESO values are indicative 


of the higher grid resolution. 


Table A.3. Grid Input. 
DUNE2D Physical Parameter 
Variable 











TFI grid generator stretching parameter 


An example of a generated grid using the transfinite interpolation technique over 


a user-defined bedform is shown in Figure (A. 1). 
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Figure A.1. . DUNE2D grid using the transfinite interpolation method. 
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4, Grid Sensitivity Study 


A grid sensitivity test of the flow performed by Andersen (1999) is summarized in 
Table (A.6). Four parameters were tested: y., M, Axcres/a, and D. The first term, y’, is 
the distance from the bottom boundary to the first grid point. Flow around the crest of 
the ripples was found to be sensitive to the bed slope. A stretched grid was applied to 
enhance the number of points in that region. Therefore, the value Axcre, the horizontal 
spacing between grid points around the crest, was tested. D is the normalized depth of 
the flow. The tests provide guidance on model performance, in particular on flow 
convergence and numerical stability. 


Table A.4. Grid Sensitivity of the Flow Test Results. 





Andersen also performed a grid sensitivity test for the suspended load. By 


varying the total vertical grid points between 40, 60, and 100, little was achieved in grid 
convergence of the suspended concentration above the bed. 


5. Model Boundary Conditions 


The lateral boundary conditions are periodic. Bottom boundary conditions are: 


a) 


=0 u,w=Oatthez=D 
Oz 
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The no-slip condition creates a viscous sub-layer. The top boundary is described by a 
symmetric boundary condition. 
6. Sediment Transport Related Setup 
_Bedload and bed concentration formula are based on Engelund and Fredsoe 
formulation, which is described in Chapter IL. Angle of repose of the sediment is 


constant at 0.65. Critical Shields parameter is 0.05 for all model runs. 

The sediment fall velocity is calculated using cline and Deigaard’s empirical 
relationship for sand using an iterative process of Equation (2.17). 

Suspended load along the lateral boundaries is treated as a true cyclic boundary. 
Additionally, wall boundary techniques are aaa to separate the suspended load from 


bedload. The bed concentration is assumed to be based on normal diffusion. 


Table A.5. | Sediment Transport Input. 


DUNE2D | Physical Parameter 





Variable 


Bedload and Bed concentration formula using Engelund-Fredsoe 
formula 
[FanL | __Non-dimensional settling velocity of he sediment ___ 
SSLO 





re Morphology Setup 


Once the morphology mode is turned on, only a few additional DUNE2D 


variables are required. The QUICK scheme with top smoothing is applied to the solution 


AS 





of the continuity equation in the morphology section below. This is third order accurate. 


(Andersen, 1999). The smoothing is based on a running average. 


Table A.6. Morphology Input. 









Physical Parameter 


[  oMscH Leonard (1979)’s Quick Scheme for continuity equation 
QSMO Smoothing is iterated 4 times over bed and at ripple crests 


B. MODEL OUTPUT 


DUNE2D 
Variable 











Figure (A.2) summarizes output parameters from each of the three modules. 







Chee 


v, LL iy 7 , a if oe Le 
WM Moye Z Uy 
iy YE, 
£4 DL f 





Computational Model 
System Equations Output 





Ys yyy YY, Ly yy , 














a LEGS LU PELE EEE “ LLIFELL LO LE. 
YY YY] WIv 
Wy YY yy) 


GEE 


ty yy Ye, 
é at + Hy YY Uy 


YY, ties Uy Yes 
yh YY ty 
Yipee yes 


UBL 34 4 04 pe “> Lis 






y, 
Wo 


ELE, 





Figure A.2. DUNE2D Module Output Diagram. 


The morphology mode simulation produces time-averaged quantities over a single 


wave period. Table (A.7) summarizes the different output variables that are saved to 
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ASCII files by DUNE2D. Most of the output is non-dimensional. Transformation of 


these quantities into meaningful dimensional output will be discussed in Appendix (C). 







* 







Table A.7. | Morphology Mode Out g 


Physical Parameter 


Bed shear,time averaged over one wave period 


Shield parameter, time averaged over one wave period 


Bedload flux, time averaged over one wave period 
Suspended load flux, time averaged over one wave period 


Total load flux, time averaged over one wave period 


Non-dimensional total load flux, time averaged over one wave 





















0 eriod 


Horizontal velocity, time averaged over one wave period 
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APPENDIX B. DUNE2D MODEL APPLICATION 


This appendix provides further information on the application of DUNE2D code. 
The first section describes general DUNE2D Fortran file information, as well as model 
‘nitiation. The following section emphasizes input file examples. 


A. FORTRAN CODE FILES 


DUNE2D is a modular Fortran 77 code program consisting of (16) *.f programs, 
(9) * programs, and a makefile. The code is divided into logical modules designed 
around the model input, flow calculations, bedload transport computations, morphology 
solutions, and output. The makefile 1s designed to compile all the programs into a single 


executable program called dune. A list of programs is given by Table (B.1). 


In the Unix environment, typing dune followed by a string four-letter project 
name, e. g. “test”, will execute the simulation application. Other files are required to run 
DUNE2D. These will be discussed in the next section. 


Table B.1. Fortran DUNE2D Files. {Program areas: a)1/O: input/output b)F: flow 
c)S: sediment transport d)M: mo G/B: grid/bounda 
du 


cmd.f: /O ft F 
if VO io.f VO 
fF mf: 

| ripple £: GB 












wave.f:F 






5] 





B. MODEL INITIATION FILES 


Once an error-free compiled version of dune exists, simulation can begin with the 
proper input files in place. The files have *.inp, *.grd, and *.fmt file extension, which 
represent input, grid, and format ASCII files. The three files must have the same four 
letter project name that is typed after dune, e. g. lest.inp, test.grd, and test.fmt. 
Additionally, special input files, Special.inp and Equations.inp, must be included in the 
simulation folder. These two files are used to initite DUNE2D simulation 


1. Input File 


The input file, est.mmp, describes the necessary input variables used to describe 
model calculation, initiation, hydrodynamic forcing, total time of integration, and other 
requirements used in flow, bedload transport, and seminions modules described in 
Chapter II. Appendix (A) provides greater detail on the specific requirements of the 
input file. An example is provided at the end of this section. 


2. Grid File 


The grid file, test.grd or test.xyb, is a two-column, ASCII file that describes the 
initial bedform grid point locations. User-defined grids are described by *.xyb files that 
are the same format as *.grd files. An example has been included. 


3. Format File 


Finally, the format file is used to describe the boundary conditions of the model. 
Since all model runs for this thesis and the work by Andersen (1999) only use a single 
block, the format file was not required in the model simulation. For more information on 


the model boundary conditions refer to Appendix (A). 
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Example of input file: test.inp 


# General commands: test.inp 03-Nov-2000 
READ, 0, 0.0 
COUT, 0, 0.0 
GourT;... Oi, 020 
SCRE, 10, 0.0 
PISO; 3; -0-0 


# Commands for the solver: 
ZTOP, 1, 0.0 

TURB, 4, 5.-87e+05 

INIT; -2,. 020 

SCHM, 4, 1.0000 

SKEW, 1, 0.0 


# Setup of the grid: 
coOoR, 1, 0.0 

GRID, 100, 40.0000 
GMET, 4, 1.2980 
RESO, 13, 0.0 


# Bed load transport: 
BEDL, 1, 0.0 

ROUG, 0, 0.0001 

TTOP, 1, 0.0863 

SPHI, 0, 0.5770 


# Suspended transport: 


SUSP, 1, 0-90 
sscy, 3, 0.0 
SSLO, 1, 0.0 
BEDC, 0, 0.0 
SMIX, 0, 0.0 
FALL, 2, 0.0168 


# Setup of the wave: 
DYNA, 1, 0.0 

ANGU, 0, 6.2832 
WAVE, 20, 1000.0000 
WEND, 0, 1.00e-24 
MEAN, 1, 0.0 


# Morphology Module: 
MORP, 1, 0.0 
MOUT, 150, 0.0 


MOIT, 5, 0.0 
MSCH, 8, 0.0 
gsmo, 4, 4.0000 
MRLX, 0, 1.0000 


# User defined variables: 
PROB, 0, 0.0 
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points. 


Partial Example of test.xyb 


120 
0.00000000 
0.00261905 
0.00523809 
0.00785714 
0.01047619 
§.01309524 
0.01571429 
Q.01833333 
0.02095238 
Q0.02357143 
0.02619048 
0.02880952 
0.03142857 
0.03404762 
0.03666667 
0.03928571 
0.04190476 
0.04452381 
0.04714286 
0.04976190 
0.05238095 
0.05500000 
0.05761905 
0.06023810 
0.06285714 
0.06547619 
0.06809524 
0.07071429 
0.07333332 
0.07595238 
0.07857143 
0.08119048 
0.08380952 
0.08642857 
0.08904762 
0.09166667 
0.09428571 
0.09690476 
0.09952381 
0.10214290 
0.10476190 
0.10738100 


40 


~0.00265241 
-0.00285512 
-0.00189552 
-0.00112030 
~0.00080800 
-0.00013576 


0.00004675 
0.00102395 
0.00214214 
0.00334799 
0.00452660 
0.00567079 
0.00685681 
0.00804007 
0.00920828 
0.01034211 
0.01155320 
0.01285784 
0.01420774 
001547332 
0.01666995 
0.01779334 
0.01893293 
0.02016495 
0.02150976 
0.02288504 
0.02416856 
0.02539719 
0.02664316 
0.02788491 
0.02913498 
0.03033831 
0.03156390 
0.03283032 
0.03408979 
0.03537232 
0.03657105 
0.03779511 
0.03902313 
0.04032772 
0.04169318 
om 





Program break. Note the first line is the total number of horizontal and vertical 


4. 


Summaries of the required files for simulation are given by Table (B. 


Summary 


files should all be placed in the same Unix folder. 
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2). These 








Table B.2. Files necessa 


Special.inp | Special input files 


Equations. inp (no user input required) 





A schematic flow chart of the input/output procedure using the required input files 


is shown in Figure (B.1). 
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Figure B.1. .DUNE2D input/output flow diagram. 


55 





THIS PAGE INTENTIALLY LEFT BLANK 





56 











APPENDIX C. DUNE2D INPUT/OUTPUT TRANSFORMATION 


This appendix provides specific input and output procedures used to implement 
DUNE2D code. An overview of the process by which data is entered into the DUNE2D 
Fortran format input files as well as transformation of the output files into meaningful, 


dimensional, as well as, graphical results is shown in Figure (C.1). 
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Figure C.1. | General input/output flow diagram of modeling simulation. 









A. GENERAL I/O PROCEDURES 


The first step to a morphology simulation focuses on determining the wave and 
current forcing dynamics. A Matlab program gathers user defined input concerning the 
ocean dynamics and sissies mean grain diameter. An *.inp file is created. Second, a 
bedform *.grd file is generated from field or experimental bedform data. A schematic 
diagram of the MATLAB programs used to generate the input file data is shown in 
Figure (C.2). 


> 


ey , Wy A f 7 
LILLE, G e ‘yj fe A fy es 

Wy lhe Ye Me EE OUI GI Ee LO LLGSEST. 
wy 


“es, Pee 
‘a % % Z x th 
GLEE AAEM AL ALE 


















Create Sub 


directories “iy, es 
Move new input fileto| 2 Le - 
This new directory Wes Le ac Se 


yey 
LLL Le 
LIAF. 


CSCI OEE 
UPLILILO CPS 


AOE fe 


(ALLL 


F 


Figure C.2. DUNE2D input data flow of MATLAB programs. 


Once all the input files are created. The first step of the morphology simulation 1s 
to mun a hot file simulation run. The main purpose of the hot file run is to generate an 
initial solution to the flow field over a flat bed. The output generated by DUNE2D 
includes a *.hot file at the end of each wave period of simulation. Normally, the hot file 
is generated over four complete wave periods. The *.hot output file is used by the actual 
morphology simulation to help initiate the model flow solutions. If the hot file is not 
included, the morphology run from a “cold” start will typically not have convergence of 


the flow solution. 


During the model simulation, output in ASCII format of specific single wave 
period averaged parameters (Table A.7) are placed into * meanbed files in the Unix 


simulation folder. MATLAB programs were developed to produce graphical data output. 


A summary of the programs and output options are included in Figure (C.2). 
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Initial bedform used as input to Dune2d 
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Figure C.2. DUNE2D output data flow of MATLAB programs. 
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les and files are provided in Figure (C.3). 
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input varia 


More detail on the specific MATLAB programs used to transform field data into 


FIELD DATA TRANSFORMATION 
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extended bedform (Figure C.4d). The filter, when applied, maintains most of the original 
roughness of the bedform except at the ends, where the smoothing takes place (Figure 
C.4e). The final step is to make the length and height of the bed non-dimensional by 
dividing through by a length scale, D (Figure C.4f). Further information on calculation 


of D can be found in Appendix (A). 
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Figure C.4. Transformation of cross-shore field bedform data to DUNE2D non- 
dimensional bedform data. 


61 














Ze Monochromatic Wave Dynamic Forcing Transformation 


MATLAB program waveforce.m is used to determine the monochromatic 
sinusoidal wave forcing described by the oscillatory velocity amplitude and peak wave 
period that are incorporated into the input file, *.imp. From the field measured u (along- 
shore) and v (cross-shore) velocity profiles of the BCDV profiler, the power spectrum can 
be found to determine the peak period, T, of the wave forcing. The oscillatory wave 
velocity amplitude is calculated as the equivalent sinusoidal velocity amplitude that gives 


the same wave variance as measured, 


(C.1) 


U, =V2Jo2 +0, 


where 0,” and o,” are the along-shore and cross-shore velocity component variances. 
This calculation is performed to get the best directional representation of the wave 
forcing over the bed. 

Finally, mean values of the w and v velocity time series represent the along-shore 


and cross-shore currents, respectively. 
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